In the Chalcogenoxide Elimination Panorama: Systematic Insight into a Key Reaction

The selenoxide elimination is a well-known reaction in organochalcogen chemistry, with wide synthetic, biological, and toxicological implications. In this work, we apply benchmarked density functional theory (DFT) calculations to investigate different aspects of the title reaction in three (bio)chemically relevant models, spanning minimal systems of theoretical interests as well as biological or synthetic organochalcogenides. The activation strain analysis (ASA) methodology is employed along a suitable reaction coordinate to obtain insight into the role of the chalcogen and of the oxidation state, to pinpoint the factors that tune the elimination reactivity of the investigated systems. Lastly, we computationally validate the hypothesis that telluroxides eliminate more slowly than selenoxides because of a detrimental hydration process that leads to unreactive hydrates.


Additional Computational Details
For the benchmark, a total of five functionals (xc), i.e., two GGAs, OLYP 1-3 and OPBE; 4 one dispersion-corrected GGA, BLYP-D3(BJ); 1-3,5-8 one hybrid, B3LYP 9,10 and one meta-hybrid M06-2X 11,12 , were preliminarily tested for the geometry optimization and energy calculations. The Slater type TZ2P basis set was used for all calculations. This basis set is of triple-ζ quality and augmented with two sets of polarization functions on each atom. For the three GGA, the small frozen core approximation was used, while for the hybrid and the metahybrid, all-electron calculations were performed since frozen core approximation is not implemented in ADF for these functionals. The role of the basis set (TZP, TZ2P and QZ4P) and of frozen core approximation (no frozen core and small core approximation) was tested for the OPBE functional, by reoptimizing all the investigated geometries and computing activation and reaction energies (Table S1) Scalar relativistic effects were included in all calculations within the zeroth-order regular approximation 13 (ZORA) as implemented in ADF. This level of theory is denoted as ZORA-xc/TZ2P(-ae). Starting from the OPBE optimized geometries (see main text), single point energies have been computed with eighteen different density functionals, i.e. ten GGAs (one dispersion-corrected GGA), two meta-GGAs, three hybrids and three meta-hybrids. In detail, BLYP, 2 BP86, 5,14 HTBS, 15 PBE, 16 mPW, 17 PW91, 18 revPBE, 19 RPBE,20 mPBE 21 were considered. In addition, the dispersion-corrected version of BP86 functional, BP86-D3(BJ), was also tested. TPSS 22,23 and SCAN 24 functionals were tested for the meta-GGAs category; PBE0, 25 OPBE0 4 and mPW1PW 17 were tested for the hybrid category (the popular B3LYP was preliminarily tested in the main text); M06, 11 M06-2X 11 and TPSSh 22 were tested for the meta-hybrid category. Frozen core (fc) approximation was not used, to allow for a rigorous comparison, since for hybrids and meta-hybrids fc is not available. All calculations are all-electron except when explicitly specified. Following this initial investigation, eighteen functionals were tested (M06-2X was included as the best performing preliminary functional, while the other four were excluded given their relatively poor performance) by running single-point energy calculations on ZORA-OPBE/TZ2P optimized geometries. All calculations were done without frozen core approximation to allow a rigorous comparison. In total, ten GGAs (one dispersion-corrected GGA), two meta-GGAs, three hybrids and three meta-hybrids were tested. The level of theory of these calculations is denoted as ZORA-xc/TZ2P-ae // ZORA-OPBE/TZ2P, and along the manuscript it will be referred to as xc // OPBE.

Extended Benchmark Discussion
The performances of DFT in reproducing CCSD(T) trends were tested as described in the additional computational details. The activation and reaction energies obtained with DFT employing the five preliminary functionals were then compared to the CCSD(T) computed reference values. The results are shown in the Table S2, while the deviation from CCSD(T) results is represented in Figure   S1. While all five functionals recover the trends discussed for CCSD(T), with the exception of the heightening of the activation energy going from Te (+2) to Te (+4) which is recovered only by OPBE and M06-2X, it can be clearly seen that the cheaper functionals (i.e. GGAs or the dispersion corrected GGA) underestimate the activation energy for the reaction of sulfoxides, selenoxides and telluroxides, with BLYP-D3(BJ) providing the worst results, with errors larger than -15 kcal mol -1 in some cases.
( Figure S1) For the other GGAs and B3LYP, the error generally increases going from S to Se, and from the lowest to the highest OS, with all the reactions in the OS +4 systematically displaying the strongest deviations from the CCSD(T) activation energies. The situation is rather different for M06-2X activation energies, that agree almost perfectly with the highly-correlated single points. With this functional, no great error arises when going from the OS 0 to the OS +4, and the ΔE ‡ of reactions involving Se shows deviation only slightly larger than those involving S.
A somewhat different picture describes deviations in reaction energies. In this case, the performance of the functionals appears to be somewhat less systematic, with some changes with the chalcogen and with the OS. Particularly, OPBE functional seems to be the worst performer for reactions involving S, but is the best performer for reactions involving Se. On the other hand, while B3LYP appears to be the worst performer for Se and Te, it is the best performer for S, with OLYP and BLYP-D3(BJ) giving similar results. In this case, M06-2X neither excels nor completely fails, predicting reaction energies within ca. ±5 kcal mol -1 with respect to CCSD(T), and always with the correct qualitative trend.
Considering these results, OLYP and OPBE functionals, benchmarked and popularly used to study SN2 reactions 4,26,27 (such as chalcogenide oxidations 28,29 ) and E2 reactions, 30,31 do not perform equally well for the quantitative description of chalcogenoxide elimination activation energies, even if they can still be used with some caution to understand the trends in the energetics in analogous elimination reactions, since CCSD(T) trends in activation and reaction energies are properly recovered also with the cheapest GGA or dispersion corrected GGA functionals.
In this preliminary analysis, OPBE appears to be the best performing GGA ( Figure S1). OPBE functional was already found to perform very well for geometry optimization of organochalcogenides Moreover, all three protocols predict activation and reaction energies that correlate very well against CCSD(T) ones, ( Figure S3) with very similar R 2 values in the range 0.97-0.99 for both activation and reaction energies and mean absolute errors of ca. 2.00 kcal mol -1 or lower for activation energies and between 2 -4 kcal mol -1 for reaction energies. Thus, in our opinion, all these three approaches can be employed to investigate the title reaction since the trends are qualitatively and quantitatively reproduced with the hybrid (OPBE0) as well as with the two meta-hybrids (M06 and M06-2X) functionals, and each of the three functionals outperforms the other two in a specific subset of reactions, with M06 being in average the best among the three.    as applied to TZ2P basis set does not seem to affect the energetics, with deviations from the TZ2Pae are only a few fractions of kcal mol -1 . Thus, TZ2P basis set is deemed to be a reasonable   compromise for the computation of chalcogenoxide elimination reactions, when small core   approximation is available.   Table S2. Activation and reaction electronic energies (kcal mol -1 ) for the for the β-elimination reaction of chalcogenoxides (minimal model) in OS 0, +2, +4. Level of theory: ZORA-xc/TZ2P (-ae)